Correlation between Counts and Clinical parameter, Microbiome: No interaction
# rnaseq clinical parameters
clinical_names <- colnames(Clin)[5:18]
clinical_names
## [1] "CD4 count (cells/ul)"
## [2] "IL-6 (pg/ml)"
## [3] "CRP (mg/ml)"
## [4] "iFABP (pg/ml)"
## [5] "sCD27 (U/ml)"
## [6] "CD14 (ng/ml)"
## [7] "LPS (pg/ml)"
## [8] "LTA\n(OD)"
## [9] "CD38+ HLA-DR+ CD4 T cells (% of CD4 T cells)"
## [10] "CD38+ HLA-DR+ CD8 T cells (% of CD8 T cells)"
## [11] "CD4 T cells (% viable)"
## [12] "CD8 T cells (% viable)"
## [13] "Tissue HIV RNA (per CD4 T cell)"
## [14] "Plasma VL"
n_clinical <- length(clinical_names)
############## rna seq ###### both
clinical_sum <- matrix(NA, 0, 7)
for (i in clinical_names[1:12]) {
output = hivornot(rnaseq, Clin, i)
# summary table
sum_a = output_sum(output, "RNAseq")
clinical_sum = rbind(clinical_sum, sum_a)
}



kable(clinical_sum, caption = "Vgenes RNAseq")
Vgenes RNAseq
| CD4 T cells (% viable) |
All Participants |
IGKV4.1 |
0.03217 |
0.0007148 |
-2.563e-02 |
32 |
| CD4 T cells (% viable) |
All Participants |
IGLV1.44 |
0.04966 |
0.0022069 |
-1.793e-02 |
32 |
| CD8 T cells (% viable) |
All Participants |
IGKV4.1 |
0.03568 |
0.0007929 |
1.901e-02 |
32 |
############## miseq ###### both
clinical_sum <- matrix(NA, 0, 7)
for (i in clinical_names[1:12]) {
output = hivornot(miseq, Clin, i)
# summary table
sum_a = output_sum(output, "MiSeq")
clinical_sum = rbind(clinical_sum, sum_a)
}










kable(clinical_sum, caption = "Vgenes MiSeq")
Vgenes MiSeq
| CD4 count (cells/ul) |
All Participants |
iga_smu |
0.03359 |
0.03359 |
4.928e-03 |
32 |
| CD4 count (cells/ul) |
All Participants |
igl_smu |
0.03720 |
0.03720 |
6.823e-03 |
32 |
| sCD27 (U/ml) |
All Participants |
iga_smu |
0.001858 |
0.001858 |
-9.638e-02 |
32 |
| CD14 (ng/ml) |
All Participants |
igl_smu |
0.01366 |
0.01366 |
-4.469e-03 |
30 |
| LPS (pg/ml) |
All Participants |
iga_smu |
0.01513 |
0.01513 |
-2.757e-01 |
30 |
| CD38+ HLA-DR+ CD4 T cells (% of CD4 T cells) |
All Participants |
iga_smu |
0.02222 |
0.02222 |
-5.492e-01 |
31 |
| CD38+ HLA-DR+ CD8 T cells (% of CD8 T cells) |
All Participants |
iga_smu |
0.0433 |
0.0433 |
-1.044e-01 |
31 |
| CD4 T cells (% viable) |
All Participants |
iga_smu |
0.001763 |
0.001763 |
1.453e-01 |
32 |
| CD8 T cells (% viable) |
All Participants |
iga_smu |
0.001005 |
0.001005 |
-1.131e-01 |
32 |
| CD8 T cells (% viable) |
All Participants |
igl_smu |
0.000425 |
0.000425 |
-1.511e-01 |
32 |
# classes of microbiome
phylum <- colnames(Clin)[1:9 + 18]
family <- colnames(Clin)[10:20 + 18]
genus <- colnames(Clin)[21:28 + 18]
species <- colnames(Clin)[29:44 + 18]
############## rna seq ###### both phylum
clinical_sum <- matrix(NA, 0, 7)
paste("phylum")
## [1] "phylum"
for (i in phylum) {
output = hivornot(rnaseq, Clin, i)
# summary table
sum_a = output_sum(output, "RNAseq")
clinical_sum = rbind(clinical_sum, sum_a)
}
kable(clinical_sum, caption = "Vgenes RNAseq phylum")
Table: Vgenes RNAseq phylum
## family
clinical_sum <- matrix(NA, 0, 7)
paste("family")
## [1] "family"
for (i in family) {
output = hivornot(rnaseq, Clin, i)
# summary table
sum_a = output_sum(output, "RNAseq")
clinical_sum = rbind(clinical_sum, sum_a)
}



kable(clinical_sum, caption = "Vgenes RNAseq family")
Vgenes RNAseq family
| Prevotellaceae |
All Participants |
IGKV1.39 |
0.01111 |
0.0002468 |
8.729e-01 |
27 |
| Brucellaceae |
All Participants |
IGLV1.40 |
0.043832 |
0.0019481 |
3.825e+02 |
27 |
| Brucellaceae |
All Participants |
IGLV3.19 |
0.007625 |
0.0001694 |
5.264e+02 |
27 |
## genus
clinical_sum <- matrix(NA, 0, 7)
paste("genus")
## [1] "genus"
for (i in genus) {
output = hivornot(rnaseq, Clin, i)
# summary table
sum_a = output_sum(output, "RNAseq")
clinical_sum = rbind(clinical_sum, sum_a)
}

kable(clinical_sum, caption = "Vgenes RNAseq genus")
Vgenes RNAseq genus
| Prevotella |
All Participants |
IGKV1.39 |
0.01966 |
0.0004368 |
8.722e-01 |
27 |
## species
clinical_sum <- matrix(NA, 0, 7)
paste("species")
## [1] "species"
for (i in species) {
output = hivornot(rnaseq, Clin, i)
# summary table
sum_a = output_sum(output, "RNAseq")
clinical_sum = rbind(clinical_sum, sum_a)
}




kable(clinical_sum, caption = "Vgenes RNAseq species")
Vgenes RNAseq species
| P. stercorea |
All Participants |
IGHV3.30 |
0.01746 |
0.0010978 |
9.537e+00 |
27 |
| P. stercorea |
All Participants |
IGKV1.39 |
0.01746 |
0.0006984 |
3.173e+00 |
27 |
| P. stercorea |
All Participants |
IGLV3.21 |
0.01746 |
0.0011640 |
5.723e+00 |
27 |
| P. oris |
All Participants |
IGHV3.73 |
0.008818 |
0.000196 |
8.854e+00 |
27 |
################## miseq #####################
clinical_sum <- matrix(NA, 0, 7)
paste("phylum")
## [1] "phylum"
for (i in phylum) {
output = hivornot(miseq, Clin, i)
# summary table
sum_a = output_sum(output, "MiSeq")
clinical_sum = rbind(clinical_sum, sum_a)
}




kable(clinical_sum, caption = "Vgenes miseq phylum")
Vgenes miseq phylum
| Cyanobacteria |
All Participants |
iga_smu |
0.00970 |
0.00970 |
1.568e+02 |
27 |
| Cyanobacteria |
All Participants |
igl_smu |
0.01856 |
0.01856 |
1.580e+02 |
27 |
| Firmicutes |
All Participants |
iga_smu |
0.008763 |
0.008763 |
8.954e+00 |
27 |
| Proteobacteria |
All Participants |
iga_smu |
0.03297 |
0.03297 |
-8.501e+00 |
27 |
## family
clinical_sum <- matrix(NA, 0, 7)
paste("family")
## [1] "family"
for (i in family) {
output = hivornot(miseq, Clin, i)
# summary table
sum_a = output_sum(output, "MiSeq")
clinical_sum = rbind(clinical_sum, sum_a)
}

kable(clinical_sum, caption = "Vgenes miseq family")
Vgenes miseq family
| Lachnospiraceae |
All Participants |
iga_smu |
0.04007 |
0.04007 |
1.007e+01 |
27 |
## genus
clinical_sum <- matrix(NA, 0, 7)
paste("genus")
## [1] "genus"
for (i in genus) {
output = hivornot(miseq, Clin, i)
# summary table
sum_a = output_sum(output, "MiSeq")
clinical_sum = rbind(clinical_sum, sum_a)
}


kable(clinical_sum, caption = "Vgenes RNAseq genus")
Vgenes RNAseq genus
| Blautia |
All Participants |
iga_smu |
0.03105 |
0.03105 |
7.152e+01 |
27 |
| Blautia |
All Participants |
igl_smu |
0.03470 |
0.03470 |
1.236e+02 |
27 |
## species
clinical_sum <- matrix(NA, 0, 7)
paste("species")
## [1] "species"
for (i in species) {
output = hivornot(miseq, Clin, i)
# summary table
sum_a = output_sum(output, "MiSeq")
clinical_sum = rbind(clinical_sum, sum_a)
}







kable(clinical_sum, caption = "Vgenes RNAseq species")
Vgenes RNAseq species
| Bacteroides_dorei |
All Participants |
iga_smu |
0.02683 |
0.02683 |
1.283e+02 |
27 |
| Bacteroides_dorei |
All Participants |
igk_smu |
0.04955 |
0.04955 |
9.830e+01 |
27 |
| Bacteroides_dorei |
All Participants |
igl_smu |
0.03370 |
0.03370 |
1.368e+02 |
27 |
| Blautia_luti |
All Participants |
iga_smu |
0.03813 |
0.03813 |
1.417e+02 |
27 |
| Blautia_luti |
All Participants |
igl_smu |
0.03847 |
0.03847 |
2.128e+02 |
27 |
| Blautia_schinkii |
All Participants |
iga_smu |
0.001614 |
0.001614 |
3.064e+02 |
27 |
| Ruminococcus_bromii |
All Participants |
iga_smu |
0.01074 |
0.01074 |
6.471e+01 |
27 |
Correlation between Counts and Clinical parameter, Microbiome: interaction with HIV
########################## linear regression ############################## equal
########################## length of outcomes and covariates
lm_Reg_Interaction <- function(gene_matrix, Clin, clin_var_name,
inter_term_name) {
######################## If there are non-numeric columns ###############
gene_matrix = gene_matrix %>% select(-HIV, -ID)
###################### Binary Interaction matrix
gene_matrix = as.matrix(gene_matrix)
# get names ready
genelistname = base::colnames(gene_matrix)
## number of gene to test, also the number of multiple test
n_gene = ncol(gene_matrix)
clinical_variable = as.matrix(Clin[, colnames(Clin) == clin_var_name])
####### the interaction variable , here a binary variable
interaction_var = as.matrix(Clin[, colnames(Clin) == inter_term_name])
ref = levels(factor(interaction_var))[1]
non_ref = levels(factor(interaction_var))[2]
interaction_var_ = factor(interaction_var, levels = c(non_ref,
ref))
## outcome lm
outcome_lm = lapply(1:n_gene, function(i) {
lm = lm(gene_matrix[, i] ~ clinical_variable * interaction_var +
Clin$Age + Clin$Gender)
## relevel to get estimate beta1 + beta3
## https://stats.stackexchange.com/questions/248248/test-whether-the-sum-of-two-coefficients-is-significantly-different-that-zero-in
lm_relevel = lm(gene_matrix[, i] ~ clinical_variable *
interaction_var_ + Clin$Age + Clin$Gender)
##### interaction term is always the last row: 6 ######## -3:
##### means no t-statistics, -2: means no Std.Error term
##### #############
coef = c(summary(lm)$coefficients[2, -c(2, 3)], summary(lm_relevel)$coefficients[2,
-c(2, 3)], summary(lm)$coefficients[6, -c(2, 3)])
return(coef)
})
outcome_lm = data.frame(matrix(unlist(outcome_lm), ncol = 6,
byrow = TRUE, dimnames = list(c(colnames(gene_matrix)),
c("slope_no", "p.value_no", "slope_hiv", "p.value_hiv",
"slope_diff", "p.value_diff"))))
# adjusted p-value
outcome_lm = outcome_lm %>% # https://cran.r-project.org/web/packages/mutoss/mutoss.pdf
# The Benjamini-Liu (BL) step-down procedure
dplyr::mutate(FDR_no = if (n_gene >= 17) {
p.adjust(p.value_no, "BH", n_gene)
} else if (n_gene <= 10) {
p.value_no
} else {
BL(p.value_no, alpha = 0.05, silent = T)$adjPValues
}, FDR_hiv = if (n_gene >= 17) {
p.adjust(p.value_hiv, "BH", n_gene)
} else if (n_gene <= 10) {
p.value_hiv
} else {
BL(p.value_hiv, alpha = 0.05, silent = T)$adjPValues
}, FDR_diff = if (n_gene >= 17) {
p.adjust(p.value_diff, "BH", n_gene)
} else if (n_gene <= 10) {
p.value_diff
} else {
BL(p.value_diff, alpha = 0.05, silent = T)$adjPValues
}, names = colnames(gene_matrix)) %>% dplyr::mutate(slope_no = round(slope_no,
10), slope_hiv = round(slope_hiv, 10), slope_diff = round(slope_diff,
10)) %>% mutate(sig = ifelse((FDR_diff < 0.05) | (FDR_hiv <
0.05), "Sig.", "Non")) %>% select(names, everything())
# sort by p.value outcome_lm =
# outcome_lm[order(outcome_lm$p.value), ]
outcome_lm = data.frame(outcome_lm)
## sample size
size = sum(!is.na(clinical_variable))
################### scatter plot ##########################
sig = sum(outcome_lm$sig == "Sig.")
if (sig > 0) {
# get significant genes list
sig_genes = outcome_lm %>% dplyr::filter(sig == "Sig.") %>%
select(names)
sig_genes = as.vector(unlist(sig_genes))
sig_genes_matrix = as.matrix(gene_matrix[, colnames(gene_matrix) %in%
sig_genes])
########## get the right sample size, missing values in clinical, and
########## miseq #########3 length(clinical_variable) =
########## max(nrow(sig_genes_matrix), length(clinical_variable) )
########## nrow(sig_genes_matrix) = max(nrow(sig_genes_matrix),
########## length(clinical_variable) )
plot_data = cbind(sig_genes_matrix, clinical_variable,
interaction_var)
colnames(plot_data) = c(sig_genes, clin_var_name, inter_term_name)
plot_data = data.frame(plot_data) %>% na.omit()
for (i in sig_genes) {
## aes
x_clin = as.numeric(plot_data[, ncol(plot_data) -
1])
y_gene = as.numeric(plot_data[, colnames(plot_data) ==
i])
plot_group = as.factor(plot_data[, ncol(plot_data)])
## with plotdata
p = ggplot(mapping = aes(x = x_clin, y = y_gene,
color = plot_group, shape = plot_group)) + geom_point() +
geom_smooth(method = "lm", se = FALSE) + labs(title = paste(sep = " ",
i, "vs.", clin_var_name, "\n group by", inter_term_name),
x = clin_var_name, y = i) + scale_color_discrete(name = paste(inter_term_name,
" Status", sep = ""), breaks = c("no", "yes"),
labels = c("Healthy Control", paste(inter_term_name,
"-infected", sep = ""))) + scale_shape_discrete(name = paste(inter_term_name,
" Status", sep = ""), breaks = c("no", "yes"),
labels = c("Healthy Control", paste(inter_term_name,
"-infected", sep = ""))) +
theme_minimal()
#
print(p)
}
} else {
}
################## summary table
return(list(results = outcome_lm, size = size, clinical = clin_var_name))
}
###################### output summary function ####################### output is
###################### the
output_sum <- function(output, genelist) {
# genelist is the name of the genelist genelist: 'RNAseq' or
# 'MiSeq' a list of two sublists
clinical = output$clinical
clinical = gsub("\n", "_", clinical, fixed = T)
# find a good file name
para_name = paste(unlist(strsplit(clinical, " ", fixed = T))[1:3],
collapse = "_")
para_name1 = gsub("/", "_", para_name, fixed = T)
size = output$size
input = output$results
# lm results
sig = sum(input$sig == "Sig.")
write.xlsx(input, paste("~/Documents/gitlab/Cario_RNASeq_Microbiom_Inte/DataProcessed/Vgenes_inter/",
sep = "", Sys.Date(), "_", para_name1, "_", genelist,
"_", "interaction.xlsx"), sheetName = paste(sep = "_",
genelist))
if (sig > 0) {
############# summary table #####################
genes = input$names[input$sig == "Sig."]
fdr_no = format.pval(input$FDR_no[input$sig == "Sig."],
digits = 4, eps = 1e-04, scientific = F)
fdr_hiv = format.pval(input$FDR_hiv[input$sig == "Sig."],
digits = 4, eps = 1e-04, scientific = F)
fdr_diff = format.pval(input$FDR_diff[input$sig == "Sig."],
digits = 4, eps = 1e-04, scientific = F)
# neg_log10p = -log10(input$p.value[input$sig == 'Sig.'])
slope_no = format(input$slope_no[input$sig == "Sig."],
digits = 4, scientific = T)
slope_hiv = format(input$slope_hiv[input$sig == "Sig."],
digits = 4, scientific = T)
slope_diff = format(input$slope_diff[input$sig == "Sig."],
digits = 4, scientific = T)
sum_a = cbind(clinical, genes, fdr_no, slope_no, fdr_hiv,
slope_hiv, fdr_diff, slope_diff, size)
colnames(sum_a) = c("Parameter", "Gene", "Adjusted_P_Non",
"Slope_Non", "Adjusted_P_HIV", "Slope_HIV", "Adjusted_P_Diff",
"Slope_Diff", "Sample_Size")
return(sum_a)
} else {
}
}
All results of RNAseq are negative.
############# with clinical viarable ################# rna seq ######
############# both
clinical_sum <- matrix(NA, 0, 9)
for (i in clinical_names[1:12]) {
output = lm_Reg_Interaction(rnaseq, Clin, i, "HIV")
# summary table
sum_a = output_sum(output, "RNAseq")
clinical_sum = rbind(clinical_sum, sum_a)
}
kable(clinical_sum, caption = "Vgenes RNAseq: Clinical Parameter and HIV Interaction")
Table: Vgenes RNAseq: Clinical Parameter and HIV Interaction
########## Miseq ################
clinical_sum <- matrix(NA, 0, 9)
for (i in clinical_names[1:12]) {
output = lm_Reg_Interaction(miseq, Clin, i, "HIV")
# summary table
sum_a = output_sum(output, "MiSeq")
clinical_sum = rbind(clinical_sum, sum_a)
}




kable(clinical_sum, caption = "Vgenes MiSeq: Clinical Parameter and HIV Interaction")
Vgenes MiSeq: Clinical Parameter and HIV Interaction
| CD4 count (cells/ul) |
igl_smu |
0.08993 |
-1.66e-02 |
0.07211 |
6.355e-03 |
0.031 |
2.295e-02 |
32 |
| sCD27 (U/ml) |
iga_smu |
0.4472 |
9.757e-02 |
0.03311 |
-7.209e-02 |
0.2038 |
-1.697e-01 |
32 |
| CD4 T cells (% viable) |
igl_smu |
0.02245 |
-2.828e-01 |
0.03657 |
4.712e-01 |
0.003171 |
7.54e-01 |
32 |
| CD8 T cells (% viable) |
igl_smu |
0.08258 |
-4.524e-01 |
0.02266 |
-1.552e-01 |
0.2616 |
2.972e-01 |
32 |
###################### with Microbiome ############ classes of microbiome
phylum <- colnames(Clin)[1:9 + 18]
family <- colnames(Clin)[10:20 + 18]
genus <- colnames(Clin)[21:28 + 18]
species <- colnames(Clin)[29:44 + 18]
####### RNA seq #################
for (j in c("phylum", "family", "genus", "species")) {
clinical_sum <- matrix(NA, 0, 9)
for (i in eval(as.name(j))) {
output = lm_Reg_Interaction(rnaseq, Clin, i, "HIV")
# summary table
sum_a = output_sum(output, "RNAseq")
clinical_sum = rbind(clinical_sum, sum_a)
}
print(kable(clinical_sum, caption = paste("Vgenes RNAseq: HIV Interaction with ",
sep = "", j)))
}






##
##
## Table: Vgenes RNAseq: HIV Interaction with phylum
##
## Parameter Gene Adjusted_P_Non Slope_Non Adjusted_P_HIV Slope_HIV Adjusted_P_Diff Slope_Diff Sample_Size
## ------------- --------- --------------- ----------- --------------- ---------- ---------------- ----------- ------------
## Fusobacteria IGKV1.39 0.9792 -4.334e+00 0.03555 9.494e+01 0.05633 9.927e+01 27
## Tenericutes IGHV4.39 0.9054 -4.584e+00 0.057447 9.071e+01 0.046978 9.530e+01 27
## Tenericutes IGKV3.20 0.7209 -1.507e+01 0.057447 1.026e+02 0.036412 1.177e+02 27
## Tenericutes IGLV2.14 0.8779 -1.095e+01 0.003595 1.829e+02 0.004242 1.939e+02 27
## Tenericutes IGLV3.19 0.9297 -1.377e+00 0.003595 1.149e+02 0.004242 1.163e+02 27
## Tenericutes IGLV6.57 0.7209 -2.309e+01 0.072953 1.427e+02 0.045772 1.658e+02 27


























##
##
## Table: Vgenes RNAseq: HIV Interaction with family
##
## Parameter Gene Adjusted_P_Non Slope_Non Adjusted_P_HIV Slope_HIV Adjusted_P_Diff Slope_Diff Sample_Size
## -------------------- --------- --------------- ----------- --------------- ----------- ---------------- ----------- ------------
## Christensenellaceae IGHV1.18 0.9416 -8.314e-01 0.0005867 8.175e+01 0.0007508 8.258e+01 27
## Christensenellaceae IGHV1.2 0.9416 -2.905e+00 0.0061953 8.499e+01 0.0047417 8.789e+01 27
## Christensenellaceae IGHV1.3 0.9416 -5.191e+00 0.0350303 9.285e+01 0.0265921 9.805e+01 27
## Christensenellaceae IGHV3.21 0.9416 1.177e+00 0.0005867 1.036e+02 0.0007508 1.025e+02 27
## Christensenellaceae IGHV3.23 0.4511 -8.260e+00 0.0072158 1.059e+02 0.0047417 1.142e+02 27
## Christensenellaceae IGHV3.7 0.9416 -4.609e+00 0.0203242 1.338e+02 0.0190335 1.385e+02 27
## Christensenellaceae IGHV3.72 0.9416 -4.629e+00 0.0515613 8.809e+01 0.0450188 9.272e+01 27
## Christensenellaceae IGHV3.74 0.5995 -7.069e+00 0.0072158 1.073e+02 0.0047417 1.144e+02 27
## Christensenellaceae IGHV4.34 0.9971 4.918e-03 0.0203242 3.688e+01 0.0205160 3.688e+01 27
## Christensenellaceae IGHV4.39 0.9416 -9.569e-01 0.0237261 5.971e+01 0.0229241 6.067e+01 27
## Christensenellaceae IGKV1.5 0.9416 -2.375e+00 0.0274026 6.646e+01 0.0240120 6.883e+01 27
## Christensenellaceae IGKV3.20 0.9416 -1.995e+00 0.0355375 6.333e+01 0.0314187 6.532e+01 27
## Christensenellaceae IGKV4.1 0.9416 -2.095e+00 0.0136861 8.262e+01 0.0122000 8.472e+01 27
## Christensenellaceae IGLV2.14 0.9971 -3.446e-01 0.0246018 9.534e+01 0.0240120 9.569e+01 27
## Christensenellaceae IGLV3.1 0.9971 3.830e-01 0.0072158 1.042e+02 0.0082207 1.039e+02 27
## Christensenellaceae IGLV3.25 0.9416 -9.937e-01 0.0216603 5.322e+01 0.0205160 5.422e+01 27
## Christensenellaceae IGLV3.9 0.9416 1.095e+00 0.0005867 1.160e+02 0.0007508 1.149e+02 27
## Christensenellaceae IGLV6.57 0.9416 -2.047e+00 0.0509046 8.875e+01 0.0450188 9.079e+01 27
## Ruminococcaceae IGLV3.9 0.713 1.385e+00 0.04973 7.179e+00 0.1447 5.794e+00 27
## Moraxellaceae IGLV6.57 0.9715 4.216e+00 0.04961 -1.695e+01 0.8657 -2.117e+01 27
## Brucellaceae IGLV3.19 0.9605 1.395e+03 0.0291 4.879e+02 0.9911 -9.074e+02 27
## Rhodospirillaceae IGHV1.18 0.9906 7.108e+00 0.009582 3.648e+03 0.009528 3.640e+03 27
## Rhodospirillaceae IGHV3.21 0.9906 1.125e+00 0.030586 4.197e+03 0.030495 4.196e+03 27
## Rhodospirillaceae IGLV2.14 0.9906 -4.838e+00 0.036897 5.092e+03 0.036449 5.096e+03 27
## Rhodospirillaceae IGLV3.25 0.9906 -2.084e+01 0.051424 2.556e+03 0.048320 2.576e+03 27
## Rhodospirillaceae IGLV3.9 0.9906 2.983e-01 0.009582 5.215e+03 0.009528 5.215e+03 27


##
##
## Table: Vgenes RNAseq: HIV Interaction with genus
##
## Parameter Gene Adjusted_P_Non Slope_Non Adjusted_P_HIV Slope_HIV Adjusted_P_Diff Slope_Diff Sample_Size
## -------------- --------- --------------- ----------- --------------- ---------- ---------------- ----------- ------------
## Coprococcus IGHV3.49 0.9165 -1.442e+01 0.0795 4.573e+01 0.03409 6.015e+01 27
## Acinetobacter IGLV6.57 0.9713 4.523e+00 0.04905 -1.7e+01 0.8436 -2.152e+01 27






##
##
## Table: Vgenes RNAseq: HIV Interaction with species
##
## Parameter Gene Adjusted_P_Non Slope_Non Adjusted_P_HIV Slope_HIV Adjusted_P_Diff Slope_Diff Sample_Size
## -------------------- --------- --------------- ---------- --------------- ----------- ---------------- ----------- ------------
## P. oris IGHV1.3 0.59822 3.718e+00 0.04368 3.594e+01 0.07205 3.223e+01 27
## P. oris IGKV3.11 0.16955 3.409e+00 0.01600 1.745e+01 0.06534 1.404e+01 27
## P. oris IGLV2.23 0.01785 6.835e+00 0.11022 -1.419e+01 0.03813 -2.103e+01 27
## P. oris IGLV3.10 0.92822 4.479e-01 0.04368 3.140e+01 0.06534 3.095e+01 27
## Blautia_glucerasei IGHV1.18 0.9941 9.382e-01 0.04646 1.168e+02 0.05263 1.158e+02 27
## Ruminococcus_bromii IGLV1.44 0.8933 -1.11e+00 0.02393 7.214e+01 0.02161 7.325e+01 27
####### miseq #################
for (j in c("phylum", "family", "genus")) {
clinical_sum <- matrix(NA, 0, 9)
for (i in eval(as.name(j))) {
output = lm_Reg_Interaction(miseq, Clin, i, "HIV")
# summary table
sum_a = output_sum(output, "MiSeq")
clinical_sum = rbind(clinical_sum, sum_a)
}
print(kable(clinical_sum, caption = paste("Vgenes MiSeq: HIV Interaction with ",
sep = "", j)))
}


##
##
## Table: Vgenes MiSeq: HIV Interaction with phylum
##
## Parameter Gene Adjusted_P_Non Slope_Non Adjusted_P_HIV Slope_HIV Adjusted_P_Diff Slope_Diff Sample_Size
## ------------ -------- --------------- ---------- --------------- ----------- ---------------- ----------- ------------
## Tenericutes iga_smu 0.04090 1.739e+02 0.01473 -6.346e+02 0.004171 -8.086e+02 27
## Tenericutes igl_smu 0.09032 1.712e+02 0.02705 -6.540e+02 0.011545 -8.252e+02 27



##
##
## Table: Vgenes MiSeq: HIV Interaction with family
##
## Parameter Gene Adjusted_P_Non Slope_Non Adjusted_P_HIV Slope_HIV Adjusted_P_Diff Slope_Diff Sample_Size
## -------------------- -------- --------------- ----------- --------------- ----------- ---------------- ----------- ------------
## Christensenellaceae iga_smu 0.5177 1.204e+01 0.01844 -4.447e+02 0.01687 -4.567e+02 27
## Rhodospirillaceae iga_smu 0.89737 2.258e+01 0.01676 -2.250e+04 0.01661 -2.252e+04 27
## Rhodospirillaceae igl_smu 0.07047 -3.773e+02 0.01766 -2.289e+04 0.01920 -2.251e+04 27




##
##
## Table: Vgenes MiSeq: HIV Interaction with genus
##
## Parameter Gene Adjusted_P_Non Slope_Non Adjusted_P_HIV Slope_HIV Adjusted_P_Diff Slope_Diff Sample_Size
## -------------- -------- --------------- ----------- --------------- ----------- ---------------- ----------- ------------
## Coprococcus iga_smu 0.04593 1.546e+02 0.006967 -3.210e+02 0.001313 -4.756e+02 27
## Coprococcus igk_smu 0.71469 2.753e+01 0.002595 -3.737e+02 0.005858 -4.013e+02 27
## Coprococcus igl_smu 0.30533 1.635e+02 0.072465 -2.777e+02 0.048302 -4.412e+02 27
## Thalassospira igl_smu 0.07232 -3.933e+02 0.03348 -2.662e+05 0.03368 -2.658e+05 27
######### 'species' ###############
clinical_sum <- matrix(NA, 0, 9)
for (i in species[-6]) {
output = lm_Reg_Interaction(miseq, Clin, i, "HIV")
# summary table
sum_a = output_sum(output, "MiSeq")
clinical_sum = rbind(clinical_sum, sum_a)
}



kable(clinical_sum, caption = paste("Vgenes MiSeq: HIV Interaction with species"))
Vgenes MiSeq: HIV Interaction with species
| Blautia_glucerasei |
iga_smu |
0.8974 |
1.364e+01 |
0.002807 |
-9.222e+02 |
0.004015 |
-9.359e+02 |
27 |
| Blautia_glucerasei |
igk_smu |
0.1062 |
-1.856e+02 |
0.033407 |
-6.528e+02 |
0.140293 |
-4.672e+02 |
27 |
| Blautia_glucerasei |
igl_smu |
0.2010 |
2.329e+02 |
0.008361 |
-9.294e+02 |
0.005640 |
-1.162e+03 |
27 |